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Abstract 

A controlled quantum system possesses a search landscape defined by the target 
physical objective as a function of the controls. This paper focuses on the landscape 
for the transition probability Pi->f between the states of a finite level quantum system. 
Traditionally, the controls are applied fields; here we extend the notion of control to also 
include the Hamiltonian structure, in the form of time independent matrix elements. 
Level sets of controls that produce the same transition probability value are shown to 
exist at the bottom -Pj->./ = 0.0 and top Pi-+f = 1.0 of the landscape with the field 
and/or Hamiltonian structure as controls. We present an algorithm to continuously 
explore these level sets starting from an initial point residing at either extreme value 
of Pi-^f- The technique can also identify control solutions that exhibit the desirable 
properties of (a) robustness at the top and (b) the ability to rapidly rise towards an 
optimal control from the bottom. Numerical simulations are presented to illustrate 
the varied control behavior at the top and bottom of the landscape for several simple 
model systems. 



1 Introduction 



Control of quantum phenomena is garnering increasing attention as growing numbers 
of experimental studies indicate the ability to find effective control over a broad range 
of quantum systems [Ij. Successful experiments using shaped laser pulses include 
selective bond breaking selective excitation of extremely similar species [3j, control 
of optical switches in semiconductors [3], selective energy transfer in bio-molecules [5], 
the manipulation of electron excitation [6], etc. 

The control landscape is defined as the physical observable as a function of the 
control variables. The general structure of the landsacpe dictates the ease of finding 
good controls as well as robustness to noise and other characteristics of a controlled 
quantum system. Thus, exploration of the landscape is important, and the tools devel- 
oped in this work are applicable to general quantum observables. For illustration, this 
paper considers the transition probability, P-i-^f, from state \i) to |/). The transition 
probability landscape has been shown [3 [H |9l [TUl E] to generally exhibit a fundamen- 
tal trap-free topology for completely controllable quantum systems, i.e., satisfaction of 
certain physical assumptions assures that there are no sub-optimal extrema capable of 
halting gradient based optimizations |12[ I13L [TS] . Large, even infinite, numbers of con- 
trol solutions are capable of achieving the same observable value defining members of 
a quantum control level set. Traditionally, the controls have taken the form of electro- 
magnetic fields applied to a fixed physical system. In contrast, a recent study treated 
the Hamiltonian structure (i.e., time independent matrix elements) as the controls 
In this circumstance, the Hamiltonian structure in the laboratory can be varied, for 
example, by substituting one functional group for another bonded to a molecular scaf- 
fold or by altering the relative composition of a material sample. Figure 1 shows the 
selection options for control variables. Depending on the control scenario, the applied 
field (Fig. la) or the Hamiltonian structure (Fig. lb) may be viewed as the control. 
The most general scenario (Fig. lc) considers both the Hamiltonian structure and the 
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applied field as dual controls. Exploiting this rich flexibility could be important in 
order to meet particularly demanding objectives (e.g., prescribed complex responses 
from a non-linear optical material). 

Control level sets, over the interior domain < -Pj_>/ < 1, have been investigated 
through a first-order (i.e., through the use of the gradient of Pj_>/ with respect to 
the controls) numerical algorithm called D-MORPH [91 1101 fT3] . The solution to the 
D-MORPH differential equation permits exploration of the connected members of the 
control level sets. The landscape top P^f = 1.0 and bottom i"?,->/ = 0.0 have partic- 
ular physical relevance. The top is the most desired location while at the bottom the 
goal is to climb away as rapidly as possible. In order to explore level sets at the top or 
bottom, this work generalizes the original first-order D-MORPH algorithm to produce 
a new second-order formulation. This extended D-MORPH algorithm utilizes the Hes- 
sian as the second derivative of Pi-^f with respect to the control variables. At both 
landscape extremes, the Hessian has a large number of zero eigenvalues and associ- 
ated eigenvectors that span the Hessian null space. The eigenvectors spanning the null 
space describe local coordinated variations of the controls that leave the Pj_>/ value 
fixed. The second-order D-MORPH algorithm systematically makes such changes of 
the controls in these null spaces to continuously explore the level sets at the landscape 
top or bottom. At this early stage of landscape feature exploration, the myopic nature 
of D-MORPH based simulations are important to link with the associated theoretical 
analyses and provide a foundation for subsequent laboratory landscape excursions with 
a variety of possible algorithms. 

The remainder of the paper is organized as follows: Section 2 presents the general 
formulation of the gradient and Hessian of the transition probability with respect to the 
controls. Utilizing the gradient and Hessian, Section 3 summarizes the algorithm used 
to explore the top and bottom of the Pi-*./ landscape. Section 4 presents numerical 
results using simple few level systems for illustration, and Section 5 offers concluding 
remarks on the significance of the work. 
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2 Transition Probability Landscape 



Throughout the paper, the symbol j3 will be used to indicate collectively the set of all 
the control parameters in the system f3j,j = 1,2, ... , M , i.e., (3 G M M . These M control 
parameters may be drawn from any component of the Hamiltonian under consideration. 
For example, these parameters may consist of the field spectral amplitudes and phases, 
as well as the Hamiltonian structure itself, i.e., the dipole matrix elements and system 
energy levels. The collective controls, f3, can be considered within the same unified 
framework regardless of their specific nature in the Hamiltonian which will be written 
as H(j3,t). In the remainder of the paper, a control that produces an interior yield 
< Pi->f < 1 will be denoted as /?, while a control producing either Pi_>/ = 0.0 or 
Pi^f = 1.0 will be written as j3* . 

The transition probability landscape is defined by Pi^f{j3) as a function of f3 



Here, U(T, 0) is the propagator evaluated at the target time T, governed by the 
Schrodinger equation 



The propagator U(t,0) is implicitly dependent on (3 through the Hamiltonian H{f3,t) 
in Eq. [2j The system has ./V states and we assume that the initial and target states, 
\i) and \f), respectively, are orthogonal, (f\i) = 0. 

Below, we construct the first and second derivatives of the propagator U(t, 0) with 
respect to the parameters /3, as they are needed for the remainder of the landscape 
analysis. On differentiating Eq. [2] with respect to an arbitrary parameter f3j, we obtain 



P^ / (/3) = |(/|[/(T,0)|i)| 2 . 



(1) 



ih^-U(t, 0) = H(J3, t)U(t, 0), U(0, 0) = 1. 



(2) 



ih-U^ (t, 0) = Hp. 03, t)U(t, 0) + H(p, t)U p . (t, 0) 



(3) 
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with initial conditions Up, (0, 0) = 0, where we use the notation Up, (t,0) = dU(t, 0) / df3j 
and Hp,((3,t) = dH(/3,t)/d/3j, . Equation [3] is an inhomogeneous equation with a 
homogeneous part identical to Eq. [2j Thus, the solution to Eq. [3J can be written as 

Up,(t,0) = (-i/H)U(t,0) f dt' U^t',O)H 0j ((3,t')U(t',O). (4) 

Jo 

The second order derivatives Up^if, 0) = d 2 U(t, 0i)/d/3id/3j can be derived similarly 
by differentiating Eq. [3] and utilizing Eqs. [2]and|H Here, we only present the result: 

U PiPj {t,0) = (-i/h)U(t,0) [ dt' tf(t',0)[H Mj (p,t')U(t',0) 

Jo 

+ H Pi {fi, t')U p , (t', 0) + Hp, (/3, t')U Pi (f , 0)], (5) 

where H^iM = 8 2 H (/3, 

2.1 Gradient of the Transition Probability Landscape 

We consider here a control /3 corresponding to a point on the landscape where Pj_^(/3) 
is not 0.0 or 1.0. The behavior of the landscape in the neighborhood of such a control 
j3 can be explored by expanding Pj_^(/3 + d/3) to first order 

P^ f {p + d/3) « P^/3) + V dPl ^ i/3) d^ (6a) 

i=i 

= P^/(/3) + VP^/(/3) T -d/3 (6b) 
where the first derivative of the transition probability, dPi->f(f3)/d/3j, is given by 

dR„ f (fl/d0 j = ^L\(f\U(T,om 2 

= (i\ul.(T,0)\f)(f\U(T,0)\i) + (i\uHT,0)\f)(WPi(T,0)\i) (?) 
= (2/h) f Im[(x\rt(t,Q)Hp.(P,t)U(t,0)\i)\dt. (8) 

JO 
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In Eq. [6bl 'T' denotes the vector transpose operation and the gradient VPj_>/(/3) will 
be used as a compact notation throughout the paper. Equation [J] was used on going 
from Eq. [7]to Eq. [8]where \\) = (T,0)\f)(f\U(T,0)\i). The first order derivative in 
Eq. [5] is useful for guiding the control j3 to either extremum Pj_>/ = 0.0 or P^f = 1.0. 

2.2 Hessian of the Transition Probability Landscape 

The top or bottom of the landscape is an extremum with a corresponding control j3* 
satisfying dPi->f(/3*)/d/3j = for j = 1, . . . , M. The behavior in the neighborhood of 
an extremum can be explored by expanding Pj^y(/3* + dfi) to second order 

-.MM „2p (o*\ 

i=l fc=i J * 

where the second derivatives of the transition probability, d 2 Pi^f{f3)/d(3idf3j, may be 
calculated by differentiating Eq. [8] with respect to /3k, leading to the general expression 

d 2 P l _ f ((3)/d(3 i dfr = 2Re{(i\uHT,0)\f)(f\U^ j (T,0)\i) 

+ Wl(T,0)\f){f\U p .(T,0)m. (10) 

Equation [10] applies at any point on the landscape and Sections 12.2.11 and 12.2.31 will, 
respectively, specialize this result to cases at the top and bottom. We can conveniently 
write Eq. [9] in matrix- vector form as 

p^fip* + dp) « p^f{n + (1/2) dp T ■ v 2 p^ f (n ■ dp. (ii) 

In Eq. [TTl V 2 Pj^/(/3*) G M. is the Hessian matrix, i.e., the matrix of second 
partial derivatives with respect to the system parameters, evaluated at f3*. 
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2.2.1 Hessian at the Top of the Landscape 

Previous formal work \17\ [T8] showed that multiple control solutions can exist at the top 
of the landscape where = 1.0. The Hessian at the top can facilitate exploration 

of the family of control solutions there. In particular, at the top of the landscape, 
U(T)\i) = e %e \f) for some phase 9 G [0, 2ir), and substituting this relationship into Eq. 
[10] upon using Eqs. 0] and [5] leads to the result 



d 2 Pj 



■ = {-■*/'<■ ) 



(-2/ft 2 )^Re[ J <p\rf(t,0)H p *(p*,t)U(t,0)\i)dt 

(12) 



o 



K (/?*)]; = ^ [ T lm[(p\uHt,0)H^(P*,t)U(t,0)\i)] dt. (13b) 



x / (i\U\t,0)Hp : (f3*,t)U(t,0)\p) dt 
Jo 

In order to reveal the structure of the Hessian, we define the set of 2N — 2 vectors: 
Vp(f3*),w p ((3*) for p = 1, . . . , i— 1, i + 1, . . . , N. Each of the vectors is of length M with 
the j-th component being 

M/3*)]j = ^ f T M(p\UHt,Q)Hp*(l3*,t)U(t,0)\i}) dt (13a) 

V2 
ft Jo 

The Hessian matrix with elements in Eq. [12] may now be written explicitly in terms of 
these vectors as 

v 2 p^ f {n = -J21 v p(p*) ■ v Kn + Mn ■ <(/?*)]• (14) 

Equation O shows that the Hessian is a sum of 2N — 2 rank one matrices. If the set 
{v p (/?*), u>p((3*)} P jti is linearly independent, then the Hessian will have rank 2N — 2; 
otherwise the rank is less than 2N—2. This result has important practical consequences 
for optimal control, as shown in the following subsections. 
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2.2.2 Robustness at the Top of the Landscape 

A physically attractive control j3* is one that consistently produces high yields even in 
the presence of some degree of noise around (3* at or near the top of the landscape. In 
the general context considered here noise may correspond to statistical variation in the 
field and/or uncertainty in the time independent Hamiltonian structure. Below, we 
show that the trace of the Hessian may be utilized as a scalar measure of the deviation 
of Pi-+f around a control /?*. 

In order to analyze the behavior at the top of the landscape, consider an arbitrary 
variation d/3 about f3* expressed as 

M 

dp = Y,db jUj (n (is) 

where {dbj}j=i,...,M are rea l expansion coefficients and the set of Uj(/3*) for j = 1, . . . M 
are the eigenvectors of V 2 Pj_ > /(/3*) satisfying 

V 2 P^/(/3*) • uj{p) = GjUjip), j = 1, . . . , M. (16) 

Importantly, only the first 2N — 2 of the eigenvalues Uj are non-zero [TTJ [18]. Here 
we assume the likely circumstance that M > 2N — 2; violation of the latter criteria will 
limit free variations on the landscape possibly leading to false traps on the otherwise 
nominally trap-free landscape. Practical circumstances will likely entail M ~ iV 2 
easily satisfying M > 2N — 2. Equation [TS] utilizes the fact that the eigenvectors of the 
Hessian form a complete set in the M-dimensional space of controls and the expansion 
coefficients in Eq. [15] may be interpreted as the projection of the control noise along 
each eigen-direction of the Hessian. The norm squared of d/3 is given by 

M 

\\dp\\ 2 ^(dp) T df3 = Y J {db 3 ) 2 . (17) 

3=1 



The effect of a control variation d/3 in Eq. [15] is given by Eq. [TT] 

M 

AP ir + f (J?) = (l/2)J2<ri (db 3 ? (18) 

2JV-2 

= (1/2) ^ a, (c%) 2 (19) 

3=1 

where we have utilized the fact that only the first 27V — 2 eigenvalues of the Hessian 
are non-zero \17\ 118] . In order to appreciate the significance of Eq. [T9l we make 
the physically reasonable assumption that the noise is equally distributed along each 
eigen-direction of the Hessian such that (dbj) 2 = (db^) 2 for all j, k = 1, . . . , M. Defining 
||d/3|| 2 = (d/3 ) 2 , then it follows from Eq. OH that dbj = ±[(d/3 ) 2 /M] 1 / 2 and Eq. Ql 
becomes 

ajw(/f) = (^jf) nv 2 p^ f (n], (20) 

since the trace of a matrix is equivalent to the sum of its eigenvalues. Utilizing Eq. [14] 
the deviation in Eq. [20] can be re- written 



p=£i j=l 



(21) 



Equations 1201 or [211 show the relationship between the transition probability deviation, 
APj_j.j(/3*), at the top of the landscape and the Hessian trace there. In the worst case 
scenario, the disturbance d/3 would lie entirely along the set of 2N — 2 eigenvectors 
corresponding to non-zero eigenvalues. In this circumstance, the factor M appearing in 
Eq. |2"01 would be replaced by 2N — 2, generally increasing the magnitude of AP^f (/?*). 
In the present work dimensionless units are used for all variables, thereby avoiding the 
need for any specific normalization of the Hessian elements. Importantly, the response 
AP^f is invariant to such a normalization. 
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2.2.3 Hessian at the Bottom of the Landscape 

Although the bottom of the landscape is not a desirable location, behavior there is 
important to explore as initial controls will likely result in -Pj->/ having a small value. 
Here, we take this to the limit and explore level sets at the absolute bottom cor- 
responding to controls (3* producing P^f = 0.0. At the bottom of the landscape, 
(i\ U^(T, 0)|/) = 0, utilizing Eqs. 4 and 5 significantly simplifies Eq. [10]to the form 

B 2 P f' T 

— ^L = (-2/h 2 )Re / (q\U^t)H^*,t)U(t)\i)dt 

x / {i\rf(t)H *(p*,t)U(t)\q) dt 
Jo 

where \q) = U<(T,Q)\f). As for the top of the landscape in Section 2.2.1, once again 
the Hessian rank at the bottom is important to assess. In order to do so, we define the 
two vectors y((3*) and z{(3*) of length M with entries 

= ^ / Re[(Q\UHt,0)Hp*(P*,t)U(t,0)\i)} dt (23a) 



o 



= ^ / lm[(q\U^t,0)Hp*(p*,t)U(t,0)\i)\ dt. (23b) 



h 

ft JO 

Then, the Hessian at the bottom may be written explicitly as 

v 2 p^ f (n = un ■ y j (n + z(n . z j ^*)] (24) 

and utilizing Eq. |24"| the Hessian trace may be expressed as 

M T 2 

Tr[V 2 P^ / (/3*)] = (2A 2 )^ / (q\uHt,0)He»(/3*,t)U(t,0)\i) dt 

3=1 Jo 



(25) 



Similar to the discussion in Section [2.2.11 regarding Eq. [T4"l now Eq. [241 shows that the 
rank of the Hessian at the bottom is at most 2. The maximal case of rank 2 occurs if 
the vectors y((3*) and z{(3*) are linearly independent; otherwise, the rank of the Hessian 
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at the bottom can be less than 2. Thus, at most there are two coordinated variations 
around the control (3* that will lift P^f(f3*) off the bottom, as discussed below. 

2.2.4 Climbing Rapidly from the Bottom of the Landscape 

At the top of the landscape Section 2.2.2 showed that the Hessian trace is an indicator 
of the robustness of a control solution. We now similarly show how the Hessian trace 
may influence gradient climbs originating from controls producing P^f values near 
the landscape bottom. To exhibit the influence of the Hessian trace on gradient based 
searches consider a control (3' near a control (3* at the bottom of the landscape. Here 
we assume that (3 1 corresponds to a point slightly off the bottom of the landscape such 
that > 0.0. The gradient at f3' = (3* + d/3 can be written as 

« V 2 ^/^*) • dfi. (26) 

In the present context, d/3 could arise either due to noise or considered as a concerted 
control variation. The gradient norm squared is 

K{p') = \\VP^ f {p')f 

^\\V 2 P^ f ((3*)-m\ 2 - (27) 

At the point (3* at the landscape bottom, the Hessian can be rewritten as 

v 2 p^ f (n = a x U1 (p*) Ul {n J + u 2 {n U2 (n T (28) 

where o\ and 01 are the two non-zero positive eigenvalues of the Hessian with U\{(3) 
and U2{(3) being the corresponding eigenvectors. Then, Eq. [271 becomes 

K((3') = a 2 (d(3 T ■ Ul ((3*)) 2 + a 2 (d(3 T ■ u 2 ((3*)) 2 . (29) 
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In order to physically interpret Eq. [29] it is convenient to consider a disturbance of 
unit magnitude \ \d(3\\ = 1. If d/3 lies entirely in the control space spanned by the two 
Hessian eigen-directions u\ and U2, then the gradient norm is bounded by a\ and 02- 
However, in the more likely scenario that the vector d/3 also has a component lying in 
the null space of the Hessian, then the gradient norm is bounded 

< ||VP^/(/3')|| < max{(Ti,ff 2 } 5 \\dj3\\ = 1. (30) 

Although maximizing the trace of the Hessian does not necessarily imply a larger upper 
bound for the gradient norm in Eq. [30| experience from numerical simulations (not 
shown) reveals that strongly curved portions of the landscape bottom generally enable 
a rapid rise in Pj_>/ when aiming at ascent. 

3 Algorithms for Landscape Explorations 

Previous formal work has shown that level sets of controls exist for all values of -Pj->/ 
including at the top and bottom of the landscape [U [TH1 [T7j . Each level set member 
produces the same transition probability value; however, controls over the level set can 
show a large variation of secondary characteristics. These characteristics may include 
the mechanism by which control is achieved, the stability of each solution to distur- 
bances, etc. Here, a D-MORPH algorithm is introduced to identify level set members 
aimed at locating control solutions displaying secondary desirable characteristics in 
addition to achieving prescribed values of -P?,->/- The myopic nature of the D-MORPH 
algorithm enables the exploration of fundamental landscape issues, as well as provides 
a link to laboratory landscape excursions. The D-MORPH method has already been 
implemented in the laboratory where it has been shown to work at the top and bottom 
of the landscape [T9J, [20]. We are interested in finding solutions f3* that belong to a 
level set denoted as C, corresponding to either Pj_^(/J*) = 1.0 or Pj^y(/3*) = 0.0. 
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In order to explore the level set, C, the control f3* is parametrized with a continuous 
parameter s, i.e., j3* /3*(s). This parameterization can be viewed as a curve f3*(s) 
as a function of s through control space. The initial control f3(s = 0) is not likely to 
lie on a level set C calling for an ascent or descent of the landscape. 

The goal of the next section is to derive differential equations for (3(s) of the form 



A particular form of Eq. [31] will produce a trajectory that can ascend/descend to 
the landscape top/bottom starting from an initial control /3 while another form of 
Eq. [31] will produce a trajectory /3*(s) G C that lies on the level set maintaining 
either Pi_>f(/3*(s)) = 1.0 or Pi^.f(f3*(s)) = 0.0 upon starting from an initial point 
at the respective landscape extremum. Section 13.11 presents the differential equation 
for ascending or descending the landscape and Section 13.21 presents the algorithm for 
exploring the level sets at values of Pj->/ = 1.0 and = 0.0. Section [331 shows 

that the D-MORPH algorithm may be utilized to locate solutions exhibiting secondary 
desirable physical characteristics along level sets at either P^t = 1.0 or 0.0. 

3.1 Ascending and Descending the Landscape 

Exploring level sets at the top or bottom of the landscape first necessitates finding an 
initial control residing at either location. For this purpose, the D-MORPH technique 
can be used to ascend or descend the landscape from an arbitrary trial control 0(0). 
Differentiating P^f(s) = Pi^f(f3(s)) with respect to s produces 



d/3(s)/ds = F[f3(s),s]. 



(31) 



dP^ f (0(s))/ds = VP^ f (f3(s)) T ■ dp(s)/ds. 



(32) 



Choosing 



dp(s)/ds = p(s)VP^ f (f3(s)) 



(33) 
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for p(s) > guarantees ascension since dPi^f(s)/ds > 0; while choosing p{s) < will 
result in descension. The function p(s) may be chosen as desired, and its magnitude 
will dictate the rate of ascent or descent. Equation [32] will be used to identify controls 
/3* residing at the landscape top and bottom. 

3.2 Movement on the Top and Bottom Level Sets 

This section is concerned with movement along the level set of control solutions either 
at the bottom or top of the control landscape, corresponding to /3*(s) values main- 
taining either Pi^.f((3*(s)) = 0.0 or Pj_ s .j(/3*(s)) = 1.0. In order to derive an equation 
for the desired trajectories, we take into account that the first derivative VPj_>/(/3*) 
is necessarily zero at the top or bottom of the landscape and therefore differentiate 
Pi^f((3*(s)) twice with respect to s, 

d 2 P^ f (fJ*(s))/ds 2 = (1/2) (d(3*(s)/ds) T ■ V 2 P^ f (P*(s)) ■ (d(3*(s)/ds) . (34) 

Investigating the trajectories /3*(s) over the top or bottom of the landscape is assured 
by demanding that d 2 P i _ ) .f((3*(s))/ds 2 = 0. 

The analysis in Sections 12.2.11 and 12.2.31 showed that the rank of the Hessian is at 
most L = 2N — 2 and L = 2 at the top and bottom, respectively. Since the Hessian 
V 2 Pi^f is a real symmetric M x M matrix, it has M real eigenvalues. The Hessian 
rank being L, with generally L < M, means that there are additionally M — L zero 
eigenvalues. Thus, the Hessian can be written as 

L 

V 2 P^ f ((3*(s)) = £>^(/3*(s)) -u, T (/?*(*)). (35) 

£=1 

where the {ue(/3* (s))} are the orthonormal eigenvectors of V 2 Pj_>./(/3*(s)) correspond- 
ing to the L non-zero eigenvalues, a^. 

Using Eqs. [341 and [35l the criterion d 2 Pi^f(s)/ds 2 = specifies the following L 
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linear constraints on d/3*(s)/ds 

u e ((3*(s)) T -dp*(s)/ds = £ = 1,...,L. (36) 

The constraints given by Eq. [36] can be re-expressed as a differential equation whose 
solution is a trajectory (3*{s) on the level set C at either the top or bottom of the 
landscape, 

d/3*(s)/ds = P c (/3*(s))-g((3*(s),s), (37) 

where g(/3*(s), s) is a vector of length M that can be freely chosen. The projector P c 
is specified by its action 

P c (f3*(s))-g(fJ*(s),s) = g(J3* (s),s)-J2 MP*(s)) lMP*(s)) T -g(P*(s),s)] (38) 

t 

which removes any component of the vector g(/3*(s),s) lying in the space spanned by 
the set of eigenvectors {ue(j3* (s))} . Thus, P c ((3*(s)) is the projector onto the null space 
of V 2 Pj_).j(/3*(s)) thereby assuring that the solution of Eq. [371 /3*(s), will satisfy the 
constraints in Eq. [36]: direct substitution of Eq. [37] into Eq. [36] upon use of Eq. [38] 
shows that the demand ug(f3* (s)) T • df3* (s) / ds = is satisfied. The choice of free vector 
g((3*(s),s) will specify the particular trajectory on the level set. 

3.3 Particular Control Trajectories at the Top or Bottom 
of the Landscape to Meet Auxiliary Objectives 

The freedom in choosing the vector g in Eq. [37] will lead to different D-MORPH 
trajectories over the top or bottom level sets. This freedom allows for secondary control 
behavior demands to be placed on level set explorations. For example, solutions may 
be sought that achieve high yield while either minimizing the effects of parameter 
disturbances on the overall outcome, minimizing the pulse fluence or optimizing other 
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physical properties. 

Consider, for example, the goal of moving on the level set C while also seeking 
to maximize a secondary cost function J(/3*), such as robustness to variations (3* — > 
P* + dfi at the top of the landscape. As the trajectory specified by the solution to Eq. 
[371 advances, changes in J will be given by 

dJ(s)/ds = VJ(/3*(s)) T • d(3*(s)/ds, (39) 

and substituting Eq. [37] into Eq. [39] gives 

dJ(s)/ds = VJ(f3*(s)) T ■ P c (f3*(s)) ■ g((3*(s),s). (40) 

As the projector P°(/3*(s)) is positive semi-definite, choosing g((3*(s), s) = p'(s)V J(/3*(s)), 
with p'(s) > will locally maximize J over the trajectory constrained to the level set 
C and the value of p'(s) will mediate the rate that J rises in the process. Similarly, 
choosing p'(s) < will lead to minimization of J over the level set C 

4 Illustrations 

Exploring the level sets at the top or bottom necessitates first either ascending or 
descending the landscape with Eq. [33] Upon reaching either P^f = 1.0 or P^f = 0.0 
(within adequate tolerance), we then move over the associated level set by solving the 
second-order D-MORPH equation in Eq. [37J The D-MORPH ascent/descent Eq. [33] 
is solved using a Runge-Kutta fourth-order variable step size integrator over the s- 
domain. The second-order D-MORPH Eq. [37] is solved using the same integrator until 
P-i^f deviates beyond a specified tolerance, 5, from the level set value. In the following 
examples, the top of the landscape is defined as Pi~^f > 1 — 5 and the bottom of the 
landscape as P%^f < d~ where 5 = 10~ 8 . If the tolerance is violated, a switch to the 
first order ascent /descent Eq. [33] is made to, respectively, ascend or descend back to 
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the level set at the top or bottom. 

The examples explore the top and bottom of the landscape for Hamiltonians of the 
form 

H(t) = H -fie(t). (41) 

where Hq is diagonal, ji is a real symmetric dipole matrix and e(t) is the electric control 
field, 

p 

e(t) = S(t) dij sm(ujijt + faj) (42) 
with P contributing frequencies and the envelope 

S(t) = ex P j- (t ~y )2 }. (43) 

The indices 'ij' in Eq. [32] correspond to the system transition frequencies and are 
explicitly stated for each particular illustration below. We set a = 3.5 in the simu- 
lations. The control variables /3 are drawn from (a) the set of Hamiltonian structure 
parameters {Ej, fj,j, /i^} where Ej is the j-th system energy level and /ij and /iy are 
diagonal and off-diagonal dipole matrix elements, respectively, and (b) the field pa- 
rameters {ctij, 4>ij} where otjj is the amplitude and (fiij the phase associated with the 
uJij spectral component. The goal is to explore controls at the top or bottom of the 
landscape through the Hamiltonian structure variables of Hq and fi in Eq. HTl and/or 
the applied field through the amplitudes and phases in Eq. [J2j In these simulations, 
dimensionless units are used, h is set to 1.0, and time is discretized between t = 
and T = 20 into 4096 points. The selected cases below are chosen to illustrate basic 
physical behavior at the top and bottom of the landscape that was identified in an 
extensive set of simulations. 
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4.1 Level Sets at the Top of the Hamiltonian Structure 
Control Landscape 

Here, we consider a simple N = 2 state system with M = 3 Hamiltonian structure 
control variables E2,fJ> 2 and fi\2 along with a fixed electric field e(t). This simple 
system is chosen for investigation as it permits a clear picture of level set exploration 
at the top of the landscape. 

For illustration, seven distinct initial systems specified by the values of -£?2) M2 and 
fj,\2 were generated with each being a point on the top of the landscape. The fixed field, 
e(t) was chosen using the form in Eq. 02] with P = 4, and e(t) merely serves the role 
of facilitating utilization of the controls E2,H2 and //i2- The field parameters are given 
in [21]. As initial systems, we chose E2 = ^2 = and a simple search over ^12 was 
performed to find that /i 12 = {0.14,0.42,0.70,0.98, 1.26, 1.54, 1.82} satisfy P^ f = 1.0 
|llj . Following the discussion in Section I2.2.1[ the rank of the Hessian at the top of 
the landscape was confirmed numerically to be 2N — 2 = 2. Thus, with three controls, 
there is an underlying one dimensional level set satisfying Pi_>2 = 1-0; the level set is a 
curve {E2(s) , (j,2(s) , fj,%2(s)} parameterized by s lying in the three dimensional control 
space. 

In the simulations, the free vector g(fi*(s),s) in Eq. [371 was chosen to be the 
eigenvector u°(/3*(s)) corresponding to the only non-zero eigenvalue of the Hessian. 
Since the eigenvectors of the Hessian are orthogonal, this choice of free vector cor- 
responds to simply moving in the direction given by u°(/3*(s)). The integration of 
Eq. [371 was performed over the domain s G [0,20]. Starting from the initial points 
{^2(0), 1^2(0), ^12(0)} above, Fig. 2 shows the seven one dimensional trajectories with 
each tracing out a path on the level set at the top of the landscape, P1-+2 = 1-0. 
Each level set trajectory is seen to be a highly coordinated path through the space 
of controls E2, \i2 and \i\2- Although it is possible that the seven paths actually are 
segments of a single overall level set trajectory (i.e., the paths eventually join together) 
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this circumstance was not found to be the case over the domain explored here. Thus, 
Fig. 2 also illustrates that level sets may have disconnected components. The highly 
contorted nature of the level set trajectories in Fig. 2 illustrates the complex structure 
of the top of the landscape. 

4.2 Robustness at the Top of the Combined Hamiltonian 
Structure and Field Control Landscape 

In this section, we explore the robustness of an N = 5 level system with fixed diagonal 
H energies E 1 = 0,E 2 = 1,E 3 = 4, E A = 9,E 5 = 16. There are M = 12 total 
controls j3* drawn from: the dipole matrix elements, {(J>j.j+i}j=i ! ...A as well as the field 
amplitudes {ojj+il^i,...^ and phases {(/>j,i+i}i=i,...,4- The remaining dipole elements 
and field amplitudes and phases were fixed at zero. The frequency associated with the 
j — > j + 1 transition is given by = |-%+i — Ej\ and the field had the form in Eq. 

02] with P = 4. 

In this illustration, the trajectory of controls /3*(s) was determined with an extra 
cost as the trace of the Hessian, 

J = Tr^P^) = f; ^ = 5> < 0, (44) 
j=i p i e=i 

introduced to explore the curvature of the level set at the top of the landscape, Pi^s = 
1.0. The value of J directly influences the robustness to variations in the controls df3 
as evident from J being the sum of the Hessian eigenvalues, An initial control at 
the top of the landscape was found through Eq. [33] with p(s) = 1.0 starting from the 
initial set of four amplitudes, phases and dipoles given in [22]. The ascent procedure 
converged to the nominal control value /3q denoting the amplitudes, phases and dipole 
values producing Pi^5 ~ 1.0. In the ascent J was not explicitly considered, but its 
value of J = —180.9 was recorded at /3q as a reference for comparison. The nominal 
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control variables, /3q, were used as the initial condition for a calculation designed to 
locate the most robust solution, /3* ob at the top of the landscape. Following the logic in 
Section [331 the free vector in Eq. [371 was chosen as g TO b(P*(s)) = VJ(/3*(s)), i.e., f3* oh 
corresponds to the smallest attained value of | J|, as J < 0. This choice of free vector 
follows from the discussion of robustness in Section [2.2.21 The gradient VJ(/3*(s)) was 
calculated using a finite difference scheme. The solution /3q, producing J = —180.9 
along with the solution /3* ob producing J = —30.3 are shown in Fig. 3a and 3b., 
respectively. 

To investigate the robustness of the control j3* at the landscape top, the Hessian 
matrix for the solutions /3q and /3* ob was calculated with Eq. [HI The Hessian eigenvalue 
of largest magnitude was —69.7 for /3q while the corresponding value was —14.1 for /3* ob . 
The eigenvector associated with the largest magnitude eigenvalue for both solutions 
is plotted in Fig 4. Fig. 4a shows that the control solution /3q is most sensitive to 
decreasing the magnitude of the dipoles /ii2 and /i23- In contrast, Fig. 4b shows that 
the control solution j3* oh is most sensitive to decreasing the magnitude of the amplitude 
045 and the dipole ^45 (the nominal value of ^45 < in Fig. 3b) as well as increasing 
the magnitude of 034 and ^34. 

Together, Figs. 3 and 4 show the manner in which a more robust solution was 
established. Over the level set trajectory (i.e., on going from Fig. 3a to Fig. 3b), the 
magnitude of the field amplitudes were decreased by ~ 70% while the magnitude of 
the dipoles were increased by approximately the same amount. This redistribution of 
parameter magnitude produces a more stable solution at the top of the landscape. The 
small magnitude of the dipole matrix elements in Fig 3a lead to very unstable solutions 
when amplitude, phase and dipole matrix element noise is introduced into the system. 
This is evident from the large eigenvector contribution of /ii2 and /i23 in Fig. 4a. 
Conversely, the balanced distribution of about equal magnitude for the field amplitudes 
and dipole elements in Fig 3b provides a more stable solution. Although the changes 
in the phases are more subtle to interpret, the optimization process manipulated these 
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as well in order to gain robustness. 

In order to assess the degree of robustness of the control solutions, a simulation 
was run to test the effect of random disturbances around /3q and /3* ob . Disturbances 
were introduced by selecting a random vector, p ran d, where each entry was chosen from 
a normal distribution with zero mean and standard deviation of 1.0. Then p ran d was 
normalized to produce ||p ran d|| = 0.25. For 5000 runs, the control solutions fifi and /3* ob 
were perturbed by the same random vector, e.g., j3' = j3* + p ran d, and Pi^{j3') was 
recorded. The results are given in Fig. 5 where panel (a) corresponds to /3q and panel 
(b) corresponds to /3* ob . 

Figure 5 shows that the nominal control solution in panel (a) is highly sensitive 
to the perturbations introduced by p ran d while the robust control solution in panel 
(b) is minimally sensitive to the random perturbations introduced by p r and- 

Although 

both control solutions are able to achieve a very high yield, ~ 1.0, only the 

solution in panel (b) can tolerate significant disturbances in the parameter settings. 
The mean value of Pi_>5 for panel (a) is 0.612 with a left standard deviation of 0.174 
while the mean value of Pi_>5 for panel (b) is 0.936 with a left standard deviation 
of 0.040. Additional tests of robustness were performed on these control solutions 
including stability to relative disturbances, i.e., (3' = j3*(l +p ran d) that confirmed the 
robustness of /3* ob over that of /3q (not shown) but to a lesser degree. 

4.3 The Ease of Climbing off the Bottom of the Control 
Field Landscape 

The bottom of the landscape is a domain where the goal is to ascend as rapidly as 
possible. As V-Pj_>j = at the bottom, the curvature dictates the rate of climbing. 
Thus, in this example, we investigate the curvature extremes encountered at the land- 
scape bottom for the target 1 — > 3 transition by analyzing the trace of the Hessian over 
families of control solutions producing Pi^3 = 0.0. An A = 4 state system is consid- 
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ered with fixed H and p specified by E\ = 0.0, E 2 = 1.0, £3 = 4.0, E^ = 9.0, p\ 2 = 
0.50, /xi3 = 0.25, /iu = 0.15,^23 = 0.20, /i24 = 0.10,^34 = 0.05, and the control field 
for variation has the form in Eq. [J2] with P = 6 corresponding to M = 12 ampli- 
tude and phase control variables. An initial control solution, /3q, was found producing 
Pi_>3 = 0.0 using the descent procedure in Eq. [33] with p(s) = —1.0, and two separate 
D-MORPH runs were initiated at /3q with the goals of either maximizing or minimizing 
the Hessian trace, J = Tr[V 2 Pj_>/], in Eq. [25] while remaining at the bottom of the 
landscape. Here, the free vector in Eq. [37] is chosen as g = ±V J, with + for maxi- 
mization and — for minimization of the Hessian trace. The control field solution that 
maximized the Hessian trace produced /S^ax with J = 8 x 10 3 , while the solution that 
minimized the Hessian trace produced /3^ in with J = 0.31. 

Figure 6 shows the dynamics of state |3) induced by both control solutions where the 
value Pi^3 — 0.0 is reached at the final time. Contrary to what might be expected, the 
population of state |3) over time is non-trivial at the landscape bottom with dynamics 
varying greatly over the breadth of the level set. In particular, P\^{t) for the maximal 
trace control solution, /3j^ ax , significantly lifts off the bottom at intermediate times to 
only finally return to — 0.0 at the final target time. In contrast, the minimal 

trace control solution, /3^ in , has P\^{t) ~ 0.0 over the full time interval < t < 20. 
Moreover, from these simulations, the mean Hessian eigenvalue at the most curved 
portion of the level set, corresponding to /3j^ ax , is found to be more than three orders 
of magnitude larger than the mean eigenvalue at the flattest portion of the level set 
corresponding to /3j^ in . 

To investigate the influence of the curvature at the landscape bottom on attempts 
to climb from there (-Pj-s-/ ~ 10~ 8 ), two gradient ascents were initiated, respectively, 
originating from the maximal and minimal trace producing fields. These were done 
over the interval s £ [0, 2] using Eq. [33] with p(s) = 1 and the results are shown in 
Fig. 7. The gradient ascent starting at the maximal trace field /^ax (corresponding 
to the left-ordinate) reaches a maximum value of Pi->3 ~ 0.999 at the final time and 
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s = 2 while the gradient ascent starting at the minimal trace field /3j^ in (corresponding 
to the right-ordinate) reaches a maximum value of P\^3 ~ 9.5 x 10~ 9 at the final time 
and s = 2. Here, and in a large number of other simulations, it is evident that the 
trace of the Hessian at the bottom of the control landscape has a significant influence 
on the rapidity with which searches can ascend from there. On the abscissa scale in 
Fig. 7, the optimization starting from the minimal trace field /3j^ in reaches a value of 
Pi^ 3 > 0.9 only at s ^ 68. 

These results on the varying degrees of curvature at the landscape bottom are 
significant because most searches for effective control fields start out with low yields. 
This circumstance is especially important for complex systems where no intuitive set 
of good trial controls may exist to begin the optimization procedure. However, these 
simulations also show that there can be attractive classes of trial controls corresponding 
to extremely curved portions of the control landscape at the bottom. 

5 Conclusion 

This work considered the transition probability landscape as a function of either the 
Hamiltonian structure and/or an applied field. Level sets were shown to exist at 
the landscape top = 1.0 and bottom Pj_>/ = 0.0. Each level set can contain 

control solutions displaying a range of secondary characteristics, including robustness 
to disturbances at the top and rapid climbing capability from the bottom. We also 
introduced an extended version of the D-MORPH algorithm for identifying desirable 
control solutions at the landscape top or bottom. This advance complements the 
prior algorithm [5] for exploring transition probability values over the interior domain 
< Pi-+f < 1.0. Taken together, these D-MORPH algorithms provide the means to 
explore the control landscape and its level sets at any value of to gain fundamental 
insights as well as for the practical concerns of finding robust control solutions. 

Although the D-MORPH procedure was developed here in a theoretical and compu- 
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tational context, it is also possible to extend some of the core concepts to the laboratory 
where statistical sampling of the controls [TO] can be used to obtain first (gradient) and 
second-order (Hessian) information (i.e., either by varying the applied field or Hamil- 
tonian structure). In particular, this procedure has been used to measure the Hessian 
at the top and bottom of the N = 4 level atomic Rubidium landscape and travel along 
the associated level sets [20] . A general procedure might implement a closed- loop tech- 
nique coupled to a learning algorithm capable of updating the field and Hamiltonian 
parameters progressively towards an optimal set at the top or bottom [23] especially 
for more complex systems. 

The ability to explore ancillary goals (e.g., robustness to field noise or Hamiltonian 
structure variations illustrated here) at the highest (and lowest) yields reveals addi- 
tional glimpses of the fundamental control landscape. Other ancillary goals can be 
envisioned including those that place demands on the nature of the mechanism that 
achieves the control goal. For example, if the target is controlled rearrangement of a 
molecule then a cost can be placed against populating any states that lead to disso- 
ciation or other undesirable physical processes. Although this study focused on the 
transition probability landscape, the concepts here can be applied to the quantum con- 
trol landscape for any observable. Furthermore, these tools can be extended to address 
an ensemble of systems using the density matrix formalism. 
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Figure 1 

Diagram of control variable options drawn from the Hamiltonian. (a) The traditional 
view of quantum control where the applied field is the control and the Hamiltonian 
structure is fixed, (b) The control is the Hamiltonian structure (e.g., energy level spac- 
ings and couplings) and the corresponding dynamics are facilitated by a fixed applied 
control field. In this case, the landscape is explored by traversing the Hamiltonian 
structure controls, (c) Both the Hamiltonian structure and the applied field are the 
controls. The overall purview in (c) reveals the richness in considering controls through 
access to molecular /material samples and external fields. With a suitable choice of con- 
trols, the entire observable landscape is subject to exploration, including the level sets 
at the top and bottom. 

Figure 2 

The P±^2 = 1 level set at the top of the Hamiltonian structure landscape for an TV = 2 
state system in Section 4.1 with controls given by by E2,^2 and /j,±2. A fixed ap- 
plied field e(t) is present facilitating the dynamics of the system. Every point along 
each trajectory corresponds to a control set E2, 1^2, M12 that accomplishes perfect yield 
under the influence of the same field. Each of the seven level set curves start with 
E2 = 0, H2 = and a particular point along the ^12 axis. 

Figure 3 

Values of the control variables for the nominal (3q and robust /3* ob control solutions for 
the example in Section 4.2 at the top of the landscape for = 1.0. These parameter 
settings correspond to field amplitudes and phases, as well as the dipole couplings as 
discussed in the text. The control /3q in (a) produces a Hessian trace of —180.9, while 
the control (3* oh in (b) produces a Hessian trace of —30.3. 
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Figure 4 

Largest eigenvector of the Hessian at the top of the landscape (Pi_>5 = 1) for the nom- 
inal and robust control solutions (a) /3g and (b) /3* ob , respectively, from the example in 
Section 4.2 and Fig. 3. In (a) the eigenvector for the nominal case has a corresponding 
eigenvalue of —69.7 while in (b) the eigenvector for the most robust case has a corre- 
sponding eigenvalue of —14.1. The nominal control solution shows sensitivity to nearly 
every dipole value. The robust control solution shows a significant sensitivity to the 
amplitudes 034 and 045, as well as the corresponding dipole elements. 

Figure 5 

Distribution of values after perturbing (a) the nominal control solution /3q and 

(b) the robust control solution /3* ob with Gaussian distributed variations of each control 
variable for the example discussed in Section 4.2 and Fig. 4. The mean Pi_>5 value 
for (a) is 0.612 with left standard deviation 0.174. Note the distinct abscissa scales 
in (a) and (b). The mean Pi_>5 value for (b) is 0.936 with left standard deviation 0.040. 

Figure 6 

Population of state |3) over time, Pi-^i), for the control field producing either a 
Hessian maximal or minimal trace at the bottom of the landscape for the example in 
Section 4.3. The minimal trace field induces dynamics that keep the population in 
state 1 3) from ever exceeding 0.03 while the maximal trace field produces much more 
complex dynamics with large population spikes over the time interval. Both control 
fields result in a final value, Pi->3 ~ 4 x 10~ 9 on the level set bottom at the final time. 

Figure 7 

Results of two separate gradient ascents of the control landscape using Eq. [33] over 
the fixed range s £ [0,2] starting near the bottom (Pi_+/ < 10~ 8 ) of the applied field 
landscape for the example in Section 4.3 and Fig. 6. The ascents start at different 
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locations on the landscape bottom corresponding to fields producing a Hessian trace 
of either maximal value (dashed: left-axis) or minimal value (solid: right-axis). Both 
initial control fields produce Pi^.3 < 1CP 8 at s = and, in both cases, the goal is to 
rapidly increase At s = 2, the ascent originating from the maximal trace field 

(trace = 8 x 10 3 ) reaches a value of Pi-^3 ~ 0.999, whereas the ascent originating from 
the minimal trace field (trace = 0.31) just reaches a value of Pi->3 ~ 9.5 x 10~ 9 . On 
the scale of the abscissa, the ascent for the minimal trace field reaches Pi-^3 > 0.9 only 
at s w 68. 
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